Synergistic use of siderophores and weak organic ligands during zinc transport in the rhizosphere controlled by pH and ion strength gradients

Citrate (Cit) and Deferoxamine B (DFOB) are two important organic ligands coexisting in soils with distinct different affinities for metal ions. It has been theorized that siderophores and weak organic ligands play a synergistic role during the transport of micronutrients in the rhizosphere, but the geochemical controls of this process remain unknown. Here we test the hypothesis that gradients in pH and ion strength regulate and enable the cooperation. To this end, first we use potentiometric titrations to identify the dominant Zn(II)–Cit and Zn(II)–DFOB complexes and to determine their ionic strength dependent stability constants between 0 and 1 mol dm−3. We parametrise the Extended Debye-Hückel (EDH) equation and determine accurate intrinsic association constants (logβ0) for the formation of the complexes present. The speciation model developed confirms the presence of [Zn(Cit)]−, [Zn(HCit)], [Zn2(Cit)2(OH)2]4−, and [Zn(Cit)2]4−, with [Zn(Cit)]− and [Zn2(Cit)2(OH)2]4− the dominant species in the pH range relevant to rhizosphere. We propose the existence of a new [Zn(Cit)(OH)3]4− complex above pH 10. We also verify the existence of two hexadentate Zn(II)–DFOB species, i.e., [Zn(DFOB)]− and [Zn(HDFOB)], and of one tetradentate species [Zn(H2DFOB)]+. Second, we identify the pH and ionic strength dependent ligand exchange points (LEP) of Zn with citrate and DFOB and the stability windows for Zn(II)–Cit and Zn(II)–DFOB complexes in NaCl and rice soil solutions. We find that the LEPs fall within the pH and ionic strength gradients expected in rhizospheres and that the stability windows for Zn(II)–citrate and Zn(II)–DFOB, i.e., low and high affinity ligands, can be distinctly set off. This suggests that pH and ion strength gradients allow for Zn(II) complexes with citrate and DFOB to dominate in different parts of the rhizosphere and this explains why mixtures of low and high affinity ligands increase leaching of micronutrients in soils. Speciation models of soil solutions using newly determined association constants demonstrate that the presence of dissolved organic matter and inorganic ligands (i.e., bicarbonate, phosphate, sulphate, or chlorides) do neither affect the position of the LEP nor the width of the stability windows significantly. In conclusion, we demonstrate that cooperative and synergistic ligand interaction between low and high affinity ligands is a valid mechanism for controlling zinc transport in the rhizosphere and possibly in other environmental reservoirs such as in the phycosphere. Multiple production of weak and strong ligands is therefore a valid strategy of plants and other soil organisms to improve access to micronutrients.

www.nature.com/scientificreports/ There is a diverse range of weak and strong binding organic ligands in the rhizosphere which originate from plants, bacteria, and fungi [1][2][3] . They are mainly secreted when there is a low supply of iron or other micronutrients in the environment 4 and they often function as allelochemicals, i.e., a ligand released by one organism directly or indirectly influences the metal homeostasis of a neighbouring organism, be it of the same species or from a different taxonomic kingdom [5][6][7] . Acting as complexing ligands, they are associated with the cycling of a wide range of biologically important micronutrients in the rhizosphere [8][9][10][11][12] . The synergistic use of strong and weak organic ligands has been scarcely explored in the microenvironment of the rhizosphere due to experimental challenges associated with the small spatial scales 13 . This knowledge gap, however, undermines our understanding of how plants, bacteria and fungi acquire micronutrients from the soil environment and take them up. This in turn hinders attempts to improve fertilizer management and biofortification in rice and other crops or to understand the reason for the multiple siderophore production observed in soil organisms 2,14 .
Ligands forming weaker complexes with metal ions include monocarboxylic, dicarboxylic, or tricarboxylic ligands. Citrate is a tricarboxylic organic ligand, which is released by rice and soil-dwelling organisms and it is present in high concentrations in soils. Ligands forming strong complexes are siderophores which are composed of a limited set of metal-binding moieties: hydroxamate, catecholate and α-hydroxycarboxylate functionalities are the most common but α-aminocarboxylate and α-hydroxyimidazole functionalities have been reported as well 15,16 . They are found at much smaller concentrations. They are typically hexadentate ligands giving octahedral coordination geometries around the central metal 17 . Desferrioxamine B (DFOB) is a hydroxamic siderophore produced by the soil bacterium Streptomyces pilosus present abundantly in soils 18 .
McRose and co-workers proposed a conceptual model for the possible interaction of different ligands released by bacteria during iron acquisition in natural solutions. This model is based upon previous mineral dissolution mechanism studies in the presence of both types of ligands 2,19,20 . The advantages of this pairing have been predicted by abiotic studies showing that a combination of siderophores and oxalate enhances Fe dissolution 19,20 . In it, as shown in Fig. 1, a ligand with low Fe affinity for example adsorbs to the mineral surface and strips the metal away from the mineral, bringing it into solution. After the formation of a labile complex at the mineral surface, a ligand exchange reaction takes place in which the metal is transferred to a ligand with high Fe affinity and the low affinity ligand is free to react with the mineral again. Modelling and experimental studies of organic ligand complexation with a wide range of first row transition metals have demonstrated that pH and ionic strength exert a critical influence on the stability and lability of metal-ligand complexes 3,[21][22][23][24] . The rhizosphere indeed is characterised by microscale gradients in pH and ionic strength, so it is reasonable to expect that these two master variables play a central part in synergistic roles of siderophores and other organic ligands during mobilisation and transport of zinc in the rhizosphere. and low Fe affinity (X) siderophore interaction 2 . Formation of a M-X complex strips M from the soil aggregate and brings the metal into solution. The metal is then exchanged between X and Y, Y transports M towards cell surface whilst X is free to return to the soil aggregate. McRose identified X as a ligand with low affinity and Y as a ligand with high affinity. To confirm the possibility of a synergistic interaction between weak and strong organic ligands with different metal binding affinities during acquisition and transport in the rhizosphere, we need to identify the physicochemical conditions in the solution that enable the exchange of metal ions from one ligand to the other, the so-called ligand exchange point (LEP), and show that distinct different stabilty windows exist. Identifying stabilty windows is also crucial to answer the uncertainty if metal-ligand complexes are stable at the relevant interfaces in the rhizosphere, i.e., soil organisms (e.g., roots, bacteria, fungi) and soil aggregates (e.g., minerals, organic matter). www.nature.com/scientificreports/ Metal complexation to ligands in solution is fast, therefore thermodynamic speciation modelling can be used to determine metal speciation and explore geochemical mechanisms that cannot be studied experimentally 25 . Geochemical speciation modelling relies on knowing the concentrations of elements in the system, the species that form (i.e., the speciation model), and accurate stability constants for the relevant complexation reactions. Thermodynamic data always refer to a selected standard state and stability constants calculated at standard state (298.15 K, 1 atm, I = 0 mol dm −3 ), where concentrations are equal to activities with no ion interactions between dissolved species, are known as intrinsic stability constants ( logβ 0 ). Stability constants measured under any other condition are conditional ( logβ ), their value depends on the chemical and physical environments (e.g., ionic strength, temperature, pressure) under which they are measured [26][27][28] . Two alternative methods are used to express the ionic medium dependence of equilibrium constants. One method takes into consideration the individual characteristics of the ionic media by using a solution dependent expression for the activity coefficient of the species involved in the equilibrium reactions. The other method uses an extended Debye Hückel expression where activity coefficients depend only on ionic charge and ionic strength, and it accounts for the solution particular properties by introducing ionic pairing between solution ions and species involved in the equilibrium reactions. The accurate determinations of logβ and logβ 0 values remain challenging and reported values can vary significantly 29,30 . For example, the stability constant reported for [Zn(Cit)]in 0.1 mol dm -3 KNO 3 differ by 1.36 log units 31,32 . A conditional stability constant can be adjusted for ionic strength and ion interactions by calculating activities of the species studied, i.e., a i = γ i [C i ], where a is the activity of the species i, γ is the activity coefficient and C is the analytically determined concentration of species i. The Davies model has been used extensively in soil studies to calculate activity coefficients of electrolytes at fairly low ionic strengths 33 . The equation has no theoretical foundation but is found to work reasonably well up to 0.1 mol dm −333 . Another widely applied approach in metal-ligand coordination studies determines conditional stability constants at a number of points within the ionic strength range of interest and uses adjustable parameters to fit the Debye-Hückel model 22,26,27 . This approach has been termed as direct method or extended Debye Hückel (EDH) model 3,22,34 . The Debye-Hückel type model is expressed by the equation where (I) is a linear function of ionic strength that can be formulated in different ways. The simplest expression for this term is f (I) = CI , with C the only adjustable parameter. This simple choice is sufficient to explain the experimental data trend generally up to 1.0 mol dm −3 .
To unravel the role of siderophores and other organic ligands during micronutrient transport in the rhizosphere is of particular interest to understand zinc efficiency mechanisms in rice genotypes 14,30,[35][36][37][38][39][40] . The pH of the rice rhizosphere can vary by more than 2 pH units within 5 mm from the root surface and tends to remain within pH 4 and 8 irrespective of the mineralogical composition of the parent material [41][42][43][44] . The approximate ionic strength of soils where rice is grown is estimated to 0.007 mol dm −345,46 , however, due to sampling difficulties, it has not yet been possible to measure ionic strength gradients directly in the rice rhizosphere. Solute concentration profiles have been simulated for rice using root exudation rates 47 . As diffusion rates can be low, simulations suggest that solute concentration increases more than 30-times within 5 mm from the root surface 47 . Assuming the ionic strength gradient parallels the solute concentration gradient (which holds if anion/cation pairs are predominantly singularly charged), ionic strength is expected to range between 0.007 and 0.21 mol dm −3 in the rice rhizosphere. When testing the possibility of cooperative interactions between high and low affinity organic ligands during zinc transport in the rhizosphere, we need to demonstrate that solution conditions exist that enable the exchange of zinc ions from one ligand to the other, so-called ligand exchange points (LEPs) 2,3 . If we find pH and/or ionic strength windows at which one Zn(II)-ligand complex becomes more stable than another and these fall within the boundaries of the gradients expected in the rhizosphere, this would provide evidence for the model of cooperative interaction as it would imply that the two types of Zn(II)-ligand complex can dominate in different parts of the rhizosphere fulfilling different functions or the same function under different soil solution condition (see Fig. 1).
The aim of the present study was to determine the pH and ion strength dependent stability fields for Zn(II) complexes with citrate (low Zn affinity) and DFOB (high Zn affinity) and the location of the ligand exchange points using accurately determined intrinsic association constants. The determination of stability constants at different ionic strengths is critical to delineate the ionic strength stability field of the ligand metal complexes 33 . The coordination chemistry of Zn(II) with citrate and DFOB is well studied from solution and solid state studies (Fig. 2). This allows us to develop an accurate speciation model and to determine relevant association constants reliably and accurately, making DFOB and citrate an ideal model system. To this end, we first study the ionic strength dependence of Zn(II)-Cit and Zn(II)-DFOB stability constants using potentiometric titrations. We set up a complete chemical speciation model for both systems, determine conditional logβ values at multiple ionic strengths in NaCl and then fit the EDH model to the data to determine the intrinsic logβ 0 using non-linear parametrisation. To our knowledge, this is the first time an accurate description of the ionic strength dependence of stability constants has been achieved for zinc in a system with a low and high affinity siderophores. We then determine the pH and ionic strength dependent LEPs for the exchange of zinc between citrate and DFOB in NaCl. We demonstrate that the predicted pH dependent LEPs fall within the gradient expected in the rice rhizosphere. Finally, we determine the pH and ionic strength dependent LEPs for the exchange of zinc between citrate and DFOB in soil solutions using the chemical composition of rice soil pore waters. We conclude that the concept of LEP and synergistic siderophore use is central to understand mobilisation and transport of zinc and other micronutrients in the rhizosphere and the reasons for multiple organic ligand production observed in soil and other environmental systems.  34 . Estimated precision was ± 0.2 mV and ± 0.003 mL for the electromotive force and titrant volume readings, respectively. All the potentiometric titrations were carried out under magnetic stirring and bubbling purified pre-saturated N 2 through the solution to exclude O 2 and CO 2 . Before studying the Zn(II)/ligand systems, the proton association constants of the ligands were determined at different ionic strengths (0.05 ≤ mol dm −3 ≤ 1.00) in NaCl solution. For each experiment, independent titrations of strong acid solutions with standard base were carried out under the same medium and ionic strength as the systems to be investigated, with the aim of determining the electrode potential (E 0 ) using the GLEE software 48 . In this way, the pH scale used was the total scale, pH = −loga H+ , where a H+ is the free proton activity. For each titration, approximately 80 to 100 data points were collected, and the equilibrium state during titrations was checked by confirming the time required to reach equilibrium.
The software programs HySS and Hyperquad were used to determine experimental conditions for the titrations and to calculate stability constants from the potentiometric data set 55,56 . For each determination of a proton association constant, hydrolysis constant or Zn-L association constant at a given ion strength and L:Zn ratio,  68 ). At pH 7 or lower, DFOB is fully protonated, including the amine, i.e., it is present as [H 4 DFOB] + ). The proposed coordination mode between Zn(II) and DFOB, based on solution studies rather than on structural data, is analogous to the one with iron, i.e., the zinc coordinates to the three hydroxamic acids, each of which is bidentate. www.nature.com/scientificreports/ at least two different titrations were performed. The titration curves for each system were treated as a single set when refining stability constants. This meant that the refinement procedure was run on both curves at the same time to derive a single set of constants. During the fitting procedure, model results with sigma values above 5 were rejected. The error reported for the conditional stability constants is the standard deviation given in the Hyperquad output file 56 . Stability constants are defined as overall association constants (β). For a polyprotic acid (H 3 L) with three acidic sites, the β constants are expressed as: Stepwise stability constants (K) are obtained from overall constants using the rule: The ionic strength model was developed by parametrising the EDH equation 57 . The logβ 0 and C parameter were determined by optimisation using the solver function from Excel following a previously developed optimisation algorithm 58 . Monte Carlo simulations were subsequently used to estimate 95% confidence intervals on the parameters of ionic strength dependence following the procedure outlined elsewhere 59 . For each species, predicted stability constants were resampled using an inverse of the cumulative normal distribution function to give sets of simulated data from which the unknown parameters were optimized. The Davies equation is used to calculate the activity coefficient using where γ i is the activity coefficient, A is a dielectric constant of the solvent, z is the charge of the ion, and I is ion strength.
The formation constant of M m L q (OH) n , logβ, determined in NaCl (1:1 salt) of the ionic strength I, is related to the corresponding value at zero ionic strength, logβ 0 by Thermodynamic and geochemical speciation modelling. Speciation calculations for the Zn(II)/citrate/DFOB systems in NaCl and in rice soil solutions were conducted using two computational speciation work packages, i.e. Hyperquad Simulation and Speciation computer software (HySS) and Visual MINTEQ 55,60,61 .
Speciation calculations with HySS used conditional and/or intrinsic stability constants determined experimentally during this study (Table 1). For all speciation calculations identifying the ligand exchange points, experiments were conducted using [Zn] = 10 -6 mol dm −3 and [L] = 10 -5 mol dm −3 , which are the standard concentrations adopted studying the effectiveness of ligands including siderophores 2,21,62 . The sensitivity of the speciation calculations was determined running computations with high and low estimates for the stability constants. The maximum error was calculated using the difference between the results of these two different speciation calculations and was 5% for speciation calculations of solutions with 1 mol dm −3 and 0.5% for speciation calculations of solutions with 0.05 mol dm −3 . We found no significance to the interpretations of the results; pH and ionic strength dependent LEPs shifted by < 0.01 pH units and < 0.001 mol dm −3 , respectively.
Speciation calculations involving solutions with the chemical composition of pore waters from rice soil in Bangladesh, including dissolved organic carbon (DOC) used Visual MINTEQ with the NICA-Donnan model to account for metal-DOC interactions 61,63,64 . The concentration of DOC was 3.3 mmol dm −3 and based on this, the concentration of fulvic acids containing carboxylic and phenolic groups are assigned by Visual MINTEQ 65 . The concentration of zinc was set at 1 × 10 -8 mol dm −366,67 . We modelled solutions with two different bicarbonate concentrations, i.e., 2 and 8 mmol dm −3 , to account for the range of alkalinity found in calcareous rice-growing soils and to study the effect of bicarbonate complexation. The citrate and DFOB concentrations tested ranged between 1 and 50 µmol dm −3 and between 0.1 and 1 µmol dm −3 , respectively. The intrinsic stability constants (logβ 0 ) for aqueous inorganic species, i.e., Zn(II)-chloride, Zn(II)-carbonate, Zn(II)-hydroxide, Zn(II)-phosphate were taken from the default database in Visual MINTEQ based on the National Institute of Standards and Technology (NIST) compilation 49 and from a comprehensive data base compilation produced by Powell and co-workers 29 . The stability constants for citrate and DFOB species determined were added to/or amended in the database. The stability constants were adjusted when running the model at each new ionic strength so that the inbuilt Davies function used in Visual MINTEQ would set logβ at the desired value, i.e., the conditional value predicted by the appropriate EDH model.
Ligand exchange points (LEPs) were calculated by determining the pH or ionic strength condition at which the concentration of Zn(II)-DFOB and Zn(II)-Cit complexes were equal.  Table 1 and SI Table 2. Also shown are the values for the optimised parameter C and the intrinsic association constants (logβ 0 ). SI Table 1 lists all the reactions included in the speciation models used to fit the potentiometric titrations and SI Fig. 2 shows single crystal X-ray structures for some of the proposed structures including ZnH 2 Cit 2 , Zn 2 Cit 2 (H 2 O) 2 and ZnCit 2 2− taken from the Cambridge Crystallographic Data Base. Figure 22 , however, a sensitivity analysis on the two Zn(II)-ligand species with the largest relative error on C found that logβ 0 remains within its error range even when logβ 0 was recalculated for the maximum and minimum possible C values. The stability constant we report for specific Zn(II)-L complexes at specific ion strengths are in line with literature reports (Fig. 3). For example, the logβ for the formation of [Zn(Cit)] − in 0.15 mol dm −3 NaCl shows good agreement with the value reported by Cigala and co-workers in 0.15 mol dm −3 NaCl; 4.79 vs. 4.71 26 . We note, however, also significant variations within reported conditional logβ values as seen Fig. 3, with published values for the formation of [Zn(HCit)] and [Zn(Cit)] − in different 1:1 electrolytes differing over two orders of magnitudes. This highlights the analytical challenges associated with accurate and precise logβ determinations of low affinity metal-ligand complexes, in low ion strength solutions 33 .
The final speciation scheme with the best statistical fits and with chemically sensible species are given in Table 1. From the eight Zn-Cit species initially considered (SI Table 1 Table 1. Ligand and proton association constants at different ionic strength (mol dm −3 , in NaCl), T = 298.1 K. The precisions for the conditional (logβ) and intrinsic (logβ 0 ) association constants are represented using experimental and model errors, respectively. The error for logβ is calculated from the titration fit using least quare optimisation within the Hyperquad program 55 . The error on the logβ 0 is calculated using Monte Carlo simulations 59 . nd not determined.  For each species, the Extended Debye-Hückel (EDH) model has been parameterised using the experimental data (see Table 1 Table 1). DFOB complexation of Zn(II) shows the same pattern of ionic strength dependence as citrate, with the greatest decrease of logβ occurring between 0.05 and 0.15 mol dm −3 NaCl, the region of most importance to the rhizosphere. The absolute decrease in [ZnL] and [Zn(HL)] stability constants between 0.05 and 0.15 mol dm −3 is approximately equal for citrate and DFOB species, average 1.58 vs. 1.73, respectively. This is explained by the effect of ionic strength primarily depending on the charge of the ions involved and free citrate and DFOB having the same electrostatic charge (−3). The ionic strength dependent parameter C shows no systematic change for neither citrate nor DFOB species. The good agreement between literature [50][51][52]54,[68][69][70] and our speciation models as well as the conditional logβ and pK a values validates the use of a single analytical method for the determination of the LEP. We note that the proposed formation of the trihydroxy Zn(II) citrate complex at pH above 10, needs to be investigated in greater detail using supplementary techniques. However, the formation of this species is not relevant for the pH range of interest in our study. As discussed below the main prevailing species in solution are those of 1:1:0 and 2:2:−2 stoichiometry for Zn:Cit:H. Figure 4 shows intrinsic stability constants for the formation of [Zn(Cit)] − and [Zn(HCit)] determined (i) using the Davies equation and the conditional association constants determined at different ionic strengths and (ii) fitting the parameterised EDH equation to the full ionic strength dataset. We find statistically significant (p < 0.05) differences for the various logβ 0 , emphasizing that the method used for activity corrections is critical when studying the small changes in speciation happening during geochemical processes in the rhizosphere 3,30,34 . Intrinsic association constants seem over predicted below 0.1 mol dm −3 and underpredicted above 0.1 mol dm −3 . The Excel file containing the calculations using Davies calculations is available on the Zenodo repository, https:// doi. org/ 10. 5281/ zenodo. 45481 62. Figure 5 shows the fraction of complexed zinc in a Zn(II)/Cit system modelled at infinite dilution using intrinsic stability constants determined (i) by fitting the EDH equation to the full citrate stability constant dataset or (ii-vi) using the Davies equation to calculate activity coefficients and adjust the citrate stability constants separately at 0.05, 0.15, 0.30, 0.50, and 1.00 mol dm −3 . At pH 5.5, the 0.05 mol dm −3 Davies-based intrinsic speciation model overpredicts the fraction of complexed zinc by approximately 20% compared to the EDH-based intrinsic ion strength model. At the same pH value, the 0.15, 0.3, 0.5, and 1 mol dm −3 Davies-based intrinsic speciation models underpredict the fraction of complexed zinc by 18, 21, 20, and 38%, respectively, compared to the EDHbased intrinsic speciation model. These results demonstrate the inconsistencies in speciation calculations that  www.nature.com/scientificreports/ arise when the same geochemical model is run using different sets of stability constants derived by applying the same method (e.g., Davies equation) to different sets of ionic strength data; even when the different sets of ionic strength data are within the activity model's ionic strength range of applicability, i.e. for Davies equation < 0.5 mol dm −3 and from the same study 3,30,34 . This is a critical observation for the pH range close to the root. A significant improvement in the accuracy of geochemical speciation calculations is achieved by adopting an empirical method with adjustable parameter when studying the ionic strength dependence of stability constants, rather than using an indirect method, as is widely practised as shown before 3,34 .

pH and ionic strength dependent ligand exchange points (LEP) for between Zn(II)-citrate and Zn(II)-DFOB and stability windows for Zn(II)-ligand complexes.
The ionic strength dependence models described above were subsequently applied to investigate the geochemical stability window of the Zn(II)/ ligand systems of interest and to determine the point of ligand exchange, i.e., where the fraction of complexed Zn(II)-Cit over fraction of complexed Zn(II)-DFOB is at unity. Figure 6 shows the fraction of complexed zinc in Zn(II)/Cit and Zn(II)/DFOB systems as a function of (a) pH and (b) ionic strength in NaCl solutions. The raw data for these plots are supplied in the supporting information (SI Table 4 and 5). We also show the range of I and pH values characteristically found between the root and soil aggregate interfaces in rice soils.

pH stability of Zn(II)-citrate and Zn(II)-DFOB complexes.
For both the Zn(II)/citrate and Zn(II)/DFOB systems, the fraction of complexed zinc increases with pH as shown in Fig. 6. For all the ionic strengths examined, Zn(II)-Cit complexes begin forming at approximately pH 3 (see SI Fig. 5 showing the species distribution for Zn(II)-Cit at different Zn:L molar ratios, ranging from 1:1 to 1:10 at 0.15 mol dm −3 NaCl solutions and different Zn and L concentrations). Once the citrate begins to bind to Zn, it takes between 6 to 7 pH units to reach total zinc complexation in the Zn(II)/Cit system. The fraction of zinc complexed by citrate increases fastest with pH in the lower ionic strength solution. At pH 6 in a 0.01 mol dm −3 www.nature.com/scientificreports/ NaCl solution, the fraction of zinc complexed by the citrate is 0.81 (SI Table 3). This is compared to just 0.07 at the same pH in a 1 mol dm −3 NaCl solution. The formation of Zn(II)-Cit complexes does not increase continuously with pH, there is a 2 to 3 pH unit plateau in the pH complexation curves for the Zn(II)/citrate system. Speciation diagrams for the Zn(II)/malate and Zn(II)/tartrate systems show similar trends to the Zn(II)/Cit pH complexation curves observed in this study; formation of Zn(II)-malate/tartrate complexes above 10% fraction of total zinc occurs at around pH 2 to 3 and there is then a plateau/only a small increase in the formation of Zn(II)-malate/tartrate complexes between pH 5 and 7 26 .
In the Zn(II)/DFOB system, ionic strength has a negligible effect on the pH complexation curves. In all solutions, complexation of zinc begins at pH 5.5 and total zinc complexation is reached within 3 pH units; no free zinc remains in the Zn(II)/DFOB systems above pH 8. The pH complexation curves for the Zn(II)/DFOB system are sigmoidal and do not contain a plateau (Fig. 6). Speciation diagrams for the Zn(II)/deoxymugineic acid (DMA), a rice produced siderophore 30 , system show a very similar pattern to the Zn(II)/DFOB pH dependent complexation curves in this study, suggesting this is a typical behaviour of high affinity siderophores. Significant concentrations of Zn(II)-DMA complexes begin forming at around pH 5 and total complexation of zinc is completed within 1.5 pH units 30 .
The pH at which Zn(II)-DFOB complexes become more abundant than Zn(II)-Cit complexes, i.e., the fraction of complexed zinc in the Zn(II)/DFOB system becomes greater than the fraction of complexed zinc in the Zn(II)/Cit system, depends strongly on ionic strength. As ionic strength increases, the pH dependent LEP becomes more acidic. In the 0.01, 0.1, and 1 mol dm −3 solutions, the pH dependent LEP is at pH 7.4, 7.1, and 6.5, respectively. This suggests that the thermodynamic favourability of the reaction for the exchange of zinc between citrate and DFOB increases with ionic strength.  Fig. 1). The raw data for these plots is supplied in the supporting information (SI Tables 3 and 4). www.nature.com/scientificreports/ For all ionic strengths tested, the predicted pH of the LEPs (pH 6.5 to 7.4) falls within the pH gradients expected in rhizospheres, which in the case of rice ranges between pH 4.0 and 8.0 40,42 . This suggests that pH gradients make it possible for Zn(II)-Cit and Zn(II)-DFOB complexes to dominate in different parts of the rice rhizosphere and, therefore, for the ligands to function synergistically. At the root interface, if pH is below 5.5 as previously suggested, neither Zn(II)-Cit nor Zn(II)-DFOB complexes are stable, excluding the possibility of an uptake of the Zn-ligand complex. It suggests that the role of the organic ligands in the uptake of metals and Zn is in increasing the leaching from soil aggregates and not in the actual uptake. Pairing of organic ligands with low and high affinity for Fe has been demonstrated by abiotic studies showing that a combination of siderophores and oxalate enhances Fe dissolution 19 .
Ionic strength stability for Zn(II)-Cit and Zn(II)-DFOB complexes. As ionic strength increases, the stability of the Zn(II)-Cit complexes decreases. Between 0 and 1 mol dm −3 ion strength, the fraction of Zn(II) complexed by citrate decreases by 0.05, 0.84, and 0.89 at pH 4, 6, and 8, respectively; in all instances this represents a relative reduction in ligand binding efficiency of approximately 92%. The Zn(II)-Cit ionic strength complexation curves initially descend sharply, two-thirds of the reduction in binding efficiency occurs before 0.2 mol dm −3 , the critical range relevant to rhizosphere gradients. For Zn(II)-DFOB complexes, ionic strength only influences the stability at pH 6. At pH 4, no Zn(II)-DFOB complexes are stable and at pH 8 Zn(II) is fully complexed by the siderophore at all ionic strength investigated. At pH 6, between 0 to 1 mol dm −3 , the decrease in fraction of complexed Zn(II) is 0.09 in the Zn(II)/DFOB system. This represents a relative reduction in ligand binding efficiency of 45%. Hence, the ionic strength controls more the stability of Zn(II)-Cit complexes than the stability of Zn(II)-DFOB complexes; it is larger (relative reduction in ligand binding efficiency 92% vs. 45%) and it is relevant over a wider pH range.
The ionic strength for the LEP at pH 6 is approximately 0.7 mol dm −3 (Fig. 6). At pH 4, the citrate remains dominant over the DFOB up to 1 mol dm −3 and at pH 8, the siderophore is already dominant over the citrate at 0 mol dm −3 . This suggests that ionic strength controlled LEP occurs at a lower ionic strength as pH increases and consequently the thermodynamic favourability of the reaction for the exchange of zinc between citrate and DFOB increases with pH. Between pH 6 and 8, the predicted LEP is at an ionic strength smaller than 0.7 mol dm −3 . Between pH 7 and 7.5, the predicted LEP ranges from ion strengths of between 0.01 and 0.1 mol dm −3 , which overlaps with the estimated ionic strength gradient expected in the rhizosphere, in the case of rice ranging between 0.01 to 0.30 mol dm −3 . Hence, our calculations suggest that when the pH of the rhizosphere is circumneutral, ionic strength gradients make it possible for Zn(II)-Cit and Zn(II)-DFOB complexes to dominate in different parts of the rhizosphere and, therefore, for the ligands to function synergistically.
Geochemical speciation calculations of the Zn/Cit/DFOB system in rice soil. Figure 7 and Table 2 show calculated LEPs for the exchange of zinc between citrate and DFOB in solutions with the chemical composition similar to that of rice soil porewater 44,45,71 . This modelling experiment aimed to assess the effect of competitive complexation by metals, dissolved organic carbon (DOC), and bicarbonate ions on the positions of the LEP and gaining insights into the importance of cooperative ligand interactions and of siderophores during zinc transport in rice soils at relevant metal and ligand concentrations as postulated previously 30,72 . The concentration of DOC used was 40 mg dm −3 . The predicted LEP from the analysis in NaCl using standard zinc and ligand (for citrate and DFOB) concentrations are included for comparison. The LEPs are reported as ranges because different sets of citrate (1 to 50 µmol dm −3 ) and DFOB (0.1 to 1 µmol dm −3 ) concentrations were analysed in solutions with 2 mmol and 8 mol dm −3 HCO 3 − , denoted as solution A and B. We find that at all ionic strengths tested, the pH of the LEPs calculated in NaCl fall within the range calculated for the pH dependent LEPs in solutions with the chemical composition similar to that of rice soils. In line with the conclusions drawn above, it is not possible to calculate the ionic strength dependent LEPs in solution A or B at pH 4 or 8. This is because for the concentration and ion strength ranges tested, citrate is dominant at pH 4 and siderophore is dominant at pH 8. The ionic strength dependent LEP predicted in NaCl solutions at pH 6 falls within the range calculated for the ionic strength dependent LEPs in the model soil solutions. The pH dependent LEPs calculated in solution A and B at individual ionic strengths vary by between 0.9 and 2 pH units depending on the concentrations of citrate and DFOB used in the speciation calculations. The ionic strength depended LEPs vary by between 0.42 and 0.44 mol dm −3 , depending on the concentrations of citrate and DFOB used in the speciation calculations. This evidence suggests that the LEPs in soil solutions are highly sensitive to ligand concentration ratios. Previous investigations have indeed highlighted the ligand concentration ratio as an important factor controlling the ligand-exchange process 73 . The size of the range calculated for the pH dependent LEPs decreases as ionic strength increases. This suggests that the ligand concentration ratio becomes less important in controlling the position of pH of the LEPs at higher ionic strength. The pH and ionic strength dependent LEPs in solution A and B are consistent i.e., the ranges broadly overlay one another. The two model solutions are differentiated by the concentration of bicarbonate ions, i.e., 2 and 8 mmol dm −3 . The agreement would imply that the LEPs are not sensitive to bicarbonate concentration. This is in line with bicarbonate ions forming significant complexes with zinc only at higher pH 29 . Indeed, numerical modelling of Zn(II) speciation with Cl -, OH -, CO 3 2-, SO 4 2-, and PO 4 3using reliable stability (formation) constants showed that in acidic and weakly alkaline freshwater systems, in the absence of organic ligands, Zn(II) speciation is dominated by Zn 2+29 . The speciation of Zn(II) is dominated by ZnCO 3 only at pH > 8. 4. In seawater systems, the speciation at pH 8.2 is dominated by Zn 2+ with ZnCl + , ZnCl 2 , Zn(CO 3 ), and Zn(SO 4 ) as minor species. Our results suggest furthermore that competition from DOC for zinc is not strong enough to shift the position of the LEPs. The role of DOC regarding possible complexation with zinc is expected to be minor given the structural properties of fulvic acids and the known low conditional stability constants for Zn(II)-DOC complexes, i.e., logK between ~ 3.8 and 4.2 74 .

Conclusions
In this study, we have developed an accurate description of the ionic strength dependence of stability constants of Zn(II)-Cit and Zn(II)-DFOB complexes and we determined the stability windows of these complexes in synthetic (NaCl) and real world (rice soil) pore water solutions. Our study allows us to draw the following important conclusions. 1   www.nature.com/scientificreports/ 2. Significant inaccuracies in speciation calculations arise when the same geochemical speciation model is run using intrinsic stability constants derived by applying the Davies equation to conditional affinity constants determined experimentally at ionic strength between 0 and 1 mol dm −3 . When the Zn(II)/Cit system is modelled at infinite dilution using intrinsic stability constants determined using adjustable parameters and the EDH model, the accuracy of geochemical speciation calculations improves by at least 18% at pH 5.5, i.e. the pH of crucial importance with respect to zinc uptake into the root.
3. For all ionic strengths examined, Zn(II)-Cit complexes begin forming at approximately pH 3. Once the citrate begins binding to zinc, it takes between 6 to 7 pH units to reach total zinc complexation in the Zn(II)/ citrate system. In the Zn(II)/DFOB system, ionic strength has a negligible effect on pH complexation curves. Complexation of zinc begins at pH 5.5 and total zinc complexation is reached within 3 pH units; no free zinc remains in the Zn(II)/DFOB systems above pH 8. 4. Variations in ionic strength affect the stability of Zn(II)-Cit complexes more than that of Zn(II)-DFOB complexes; it is larger (relative reduction in ligand binding efficiency is twofold) and it is relevant over a wider pH range. The pH range at which Zn(II)-DFOB complexes become more stable than Zn(II)-citrate complexes depends critically on the ionic strength. Ion strength gradients within the circumneutral pH of the rhizosphere make it possible for Zn(II)-Cit and Zn(II)-DFOB complexes to dominate in different parts of the rhizosphere and, therefore, for the ligands to function synergistically. This likely explains the previously observed synergistic effect of siderophores and weak organic acids enhancing Fe dissolution.
5. pH and ionic strength dependent LEPs for citrate and DFOB with zinc in NaCl solutions and in solutions based on the chemistry of rice soils fall within the pH and ionic strength gradients expected in the rhizosphere of rice soils. The LEPs in soil solutions are highly sensitive to ligand concentration ratios, in line with the idea that ligand concentration ratios play an control of ligand-exchange processes. The presence of DOC and inorganic ligands has no significant effect.
6. Synergistic use of weak and strong binding organic ligands is an important mechanism for rice and other soil organisms to overcome micronutrient deficiencies. The synergistic use of organic ligands with high and low Zn affinity is possible due to pH and ion strength gradients found in the rhizospheres. Siderophores can have a distinct non-classical function to take up metals other than Fe and play an important role as zincophores in rice soils. In line with recent findings 3,24,30,34,36,43,72 , the present study supports the conclusion that siderophores are central to access Zn in the soil.